// ==================================================================================
// Shared Genomics Project MPI Codebase
// Version 1.0 30/04/2010
//
// (c) 2010 University of Manchester all rights reserved
//
// This file is distributed under the GNU General Public License, Version 2.  
// Please see the file COPYING.txt for more details
// ==================================================================================

/*!
\file
\brief Implementation of the SNP Analysis Examples 
*/

#include <stdio.h>
#include <stdlib.h>
#include <CUnit.h>

#include "examples.h"
#include "commands.h"
#include "SMG_All.h"

void echo_spacer(void) {
	printf("\n========================================\n");
}

void echo_command(char* cmd) {
	printf("\nExecuting: %s\n", cmd);
}

void ex_0001_binary_encoding_dataset(void) {
	int result;
	char command[1000];
	char *map = "../data/snp_1000.map"; // Map file (SNP list)
	char *ped = "../data/snp_1000.ped"; // Pedigree File (Allelic Data)
	char *id = "ex_0001"; // binary encoded data-set ID

	// Create the source data-set.
	// This is a one-off operation with the binary encoded data-set persisting somewhere on 
	// the file system.
	fill_binary_ped_command(command, map, ped, id, LDIR, RDIR);
	echo_command(command);
	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	printf("\nFile Outputs:-\n");	
	printf("\nex_0001.bmap.dat = Binary Encoded Map File (SNP List)\n");
	printf("ex_0001.bped.dat = Binary PED File (Boolean Data)\n");
	printf("ex_0001.ff.dat = Founders Data (gender information)\n");
	printf("ex_0001.pm.dat = Person Missingness for each individual in Data-set\n");
	printf("ex_0001.summary.dat = Summary statistics of SNP/Person missingness rates\n");
	echo_spacer();
}

void ex_0002_filtering_binary_encoding_dataset_s(void) {
	int result;
	char command[1000];

	// Text PLINK-compatible input files.
	char *map = "../data/snp_1000.map"; // Map file (SNP list)
	char *ped = "../data/snp_1000.ped"; // Pedigree File (Allelic Data)
	char *phe = "../data/random_phe_data.txt"; // Phenotype data

	// Filter paramaters
	char *rowMajor = "s"; // SNPs
	double mind = 0.03; // Missingness threshold person
	double geno = 0.031; // SNP missingness rate
	double maf = 0.15; // Major allelic frequency
	double hwe = 0.10; // HWE inclusion threshold
	int mode = 1; // Filtering Step (1 or 2)
	int nperm = 10000; // No. permutations required

	// Files of the binary encoded data-set created by "binary_ped.exe"
	char *id = "ex_0002"; // binary encoded data-set ID
	char *run_id = "run_ex_0002"; // Filtered data-set ID
	char *bmap = "ex_0002.bmap.dat"; // Data-set SNP list
	char *bped = "ex_0002.bped.dat"; // Data-set allelic data
	char *person_missingness = "ex_0002.pm.dat"; // Person missingness data for unfiltered data-set
	char *founders = "ex_0002.ff.dat";  // Gender Data
	char *phenoName = "BT";

	// Create the source data-set.
	// This is a one-off operation with the binary encoded data-set persisting somewhere on 
	// the file system.
	fill_binary_ped_command(command, map, ped, id, LDIR, RDIR);
	echo_command(command);
	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// First call of filter program to process and filter an allelic data-set.
	mode = 1;
	fill_pformat_data2_command(command, NUMPROCS, bmap, bped, phe, phenoName, 
							   NUMPROCS, run_id, LDIR, RDIR, rowMajor, nperm, mode, mind, geno, maf,
							   person_missingness, founders, hwe);

	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// 2nd call to join the fragments and create the filtered data-set.
	// This call is always a single core operation.
	mode = 2;
	fill_pformat_data2_command(command, 1, bmap, bped, phe, phenoName, 
							   NUMPROCS, run_id, LDIR, RDIR, rowMajor, nperm, mode, mind, geno, maf,
							   person_missingness, founders, hwe);

	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// On a 560K SNP data-set of 1000, 32 seconds wre required to loop through the 
	// data-set and a further 10 seconds to write the output file.
	// Running a filter job on 4 cores + improvements to the file I/O
	// dropped the wait time to 15 seconds for a filter job to complete.

	printf("\nFile Outputs:-\n");
	printf("\nrun_ex_0002.bmap.dat = Filtered SNP List\n");
	printf("run_ex_0002.fslog.txt = log file containing filter settings\n");
	printf("run_ex_0002.perm.dat = Permutation indices file\n");
	printf("run_ex_0002.ptd.dat = Persons to delete list for strata generation\n");
	printf("run_ex_0002.rseed.dat = Array of Random seeds for cores on a cluster\n");
	printf("run_ex_0002.rseed.txt = Initialisation value used by PRNG\n");
	printf("run_ex_0002.samples.txt = Sample array containg gender and phenotype data\n");
	printf("run_ex_0002.sped.dat = BitString array of allelic data (SNP major)\n");	
	echo_spacer();
}

void ex_0003_filtering_binary_encoding_dataset_i(void) {
	int result;
	char command[1000];

	// Text PLINK-compatible input files.
	char *map = "../data/snp_1000.map"; // Map file (SNP list)
	char *ped = "../data/snp_1000.ped"; // Pedigree File (Allelic Data)
	char *phe = "../data/random_phe_data.txt"; // Phenotype data

	// Filter paramaters
	char *rowMajor = "i"; // Individuals
	double mind = 0.03; // Missingness threshold person
	double geno = 0.031; // SNP missingness rate
	double maf = 0.15; // Major allelic frequency
	double hwe = 0.10; // HWE inclusion threshold
	int mode = 1; // Filtering Step (1 or 2)
	int nperm = 10000; // No. permutations required

	// Files of the binary encoded data-set created by "binary_ped.exe"
	char *id = "ex_0003"; // binary encoded data-set ID
	char *run_id = "run_ex_0003"; // Filtered data-set ID
	char *bmap = "ex_0003.bmap.dat"; // Data-set SNP list
	char *bped = "ex_0003.bped.dat"; // Data-set allelic data
	char *person_missingness = "ex_0003.pm.dat"; // Person missingness data for unfiltered data-set
	char *founders = "ex_0003.ff.dat"; // Gender Data
	char *phenoName = "BT"; // Column to read from the phenotype file

	// Create the source data-set.
	// This is a one-off operation with the binary encoded data-set persisting somewhere on 
	// the file system.
	fill_binary_ped_command(command, map, ped, id, LDIR, RDIR);
	echo_command(command);
	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// First call of filter program to process and filter an allelic data-set.
	mode = 1;
	fill_pformat_data2_command(command, NUMPROCS, bmap, bped, phe, phenoName, 
							   NUMPROCS, run_id, LDIR, RDIR, rowMajor, nperm, mode, mind, geno, maf,
							   person_missingness, founders, hwe);

	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// 2nd call to join the fragments and create the filtered data-set.
	// This call is always a single core operation.
	mode = 2;
	fill_pformat_data2_command(command, 1, bmap, bped, phe, phenoName, 
							   NUMPROCS, run_id, LDIR, RDIR, rowMajor, nperm, mode, mind, geno, maf,
							   person_missingness, founders, hwe);

	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// On a 560K SNP data-set of 1000, 32 seconds wre required to loop through the 
	// data-set and a further 10 seconds to write the output file.
	// Running a filter job on 4 cores + improvements to the file I/O
	// dropped the wait time to 15 seconds for a filter job to complete.

	printf("\nFile Outputs:-\n");
	printf("\nrun_ex_0002.bmap.dat = Filtered SNP List\n");
	printf("run_ex_0002.fslog.txt = log file containing filter settings\n");
	printf("run_ex_0002.perm.dat = Permutation indices file\n");
	printf("run_ex_0002.ptd.dat = Persons to delete list for strata generation\n");
	printf("run_ex_0002.rseed.dat = Array of Random seeds for cores on a cluster\n");
	printf("run_ex_0002.rseed.txt = Initialisation value used by PRNG\n");
	printf("run_ex_0002.samples.txt = Sample array containg gender and phenotype data\n");
	printf("run_ex_0002.iped.dat = BitString array of allelic data (Individual major)\n");

	printf("\nN/B Individual major requires a rotation of the allelic PED array.\n");
	printf("Does a straight x<->y swap by coping between pre-allocated arrays.\n");
	printf("Do not use this function on massive data-sets.\n");
	echo_spacer();
}

void ex_0004_filtering_binary_encoding_dataset_s(void) {
	int result;
	char command[1000];

	// Text PLINK-compatible input files.
	char *map = "../data/snp_1000.map"; // Map file (SNP list)
	char *ped = "../data/snp_1000.ped"; // Pedigree File (Allelic Data)
	char *phe = "../data/random_phe_data.txt"; // Phenotype data

	// Filter paramaters
	char *rowMajor = "s"; // SNPs
	double mind = 0.03; // Missingness threshold person
	double geno = 0.031; // SNP missingness rate
	double maf = 0.15; // Major allelic frequency
	double hwe = 0.10; // HWE inclusion threshold
	int mode = 1; // Filtering Step (1 or 2)
	int nperm = 10000; // No. permutations required

	// Files of the binary encoded data-set created by "binary_ped.exe"
	char *id = "ex_0004"; // binary encoded data-set ID
	char *run_id = "run_ex_0004"; // Filtered data-set ID
	char *bmap = "ex_0004.bmap.dat"; // Data-set SNP list
	char *bped = "ex_0004.bped.dat"; // Data-set allelic data
	char *person_missingness = "ex_0004.pm.dat"; // Person missingness data for unfiltered data-set
	char *founders = "ex_0004.ff.dat";	// Gender Data
	char *phenoName = "QT"; // Column to read from the phenotype file

	// Create the source data-set.
	// This is a one-off operation with the binary encoded data-set persisting somewhere on 
	// the file system.
	fill_binary_ped_command(command, map, ped, id, LDIR, RDIR);
	echo_command(command);
	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// First call of filter program to process and filter an allelic data-set.
	mode = 1;
	fill_pformat_data2_command(command, NUMPROCS, bmap, bped, phe, phenoName, 
							   NUMPROCS, run_id, LDIR, RDIR, rowMajor, nperm, mode, mind, geno, maf,
							   person_missingness, founders, hwe);

	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// 2nd call to join the fragments and create the filtered data-set.
	// This call is always a single core operation.
	mode = 2;
	fill_pformat_data2_command(command, 1, bmap, bped, phe, phenoName, 
							   NUMPROCS, run_id, LDIR, RDIR, rowMajor, nperm, mode, mind, geno, maf,
							   person_missingness, founders, hwe);

	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// On a 560K SNP data-set of 1000, 32 seconds wre required to loop through the 
	// data-set and a further 10 seconds to write the output file.
	// Running a filter job on 4 cores + improvements to the file I/O
	// dropped the wait time to 15 seconds for a filter job to complete.

	printf("\nFile Outputs:-\n");
	printf("\nrun_ex_0004.bmap.dat = Filtered SNP List\n");
	printf("run_ex_0004.fslog.txt = log file containing filter settings\n");
	printf("run_ex_0004.perm.dat = Permutation indices file\n");
	printf("run_ex_0004.ptd.dat = Persons to delete list for strata generation\n");
	printf("run_ex_0004.rseed.dat = Array of Random seeds for cores on a cluster\n");
	printf("run_ex_0004.rseed.txt = Initialisation value used by PRNG\n");
	printf("run_ex_0004.samples.txt = Sample array containg gender and phenotype data\n");
	printf("run_ex_0004.sped.dat = BitString array of allelic data (SNP major)\n");
	echo_spacer(); 
}

void ex_0005_filtering_binary_encoding_dataset_i(void) {
	int result;
	char command[1000];

	// Text PLINK-compatible input files.
	char *map = "../data/snp_1000.map"; // Map file (SNP list)
	char *ped = "../data/snp_1000.ped"; // Pedigree File (Allelic Data)
	char *phe = "../data/random_phe_data.txt"; // Phenotype data

	// Filter paramaters
	char *rowMajor = "i"; // Individuals
	double mind = 0.03; // Missingness threshold person
	double geno = 0.031; // SNP missingness rate
	double maf = 0.15; // Major allelic frequency
	double hwe = 0.10; // HWE inclusion threshold
	int mode = 1; // Filtering Step (1 or 2)
	int nperm = 10000; // No. permutations required

	// Files of the binary encoded data-set created by "binary_ped.exe"
	char *id = "ex_0005"; // binary encoded data-set ID
	char *run_id = "run_ex_0005"; // Filtered data-set ID
	char *bmap = "ex_0005.bmap.dat"; // Data-set SNP list
	char *bped = "ex_0005.bped.dat"; // Data-set allelic data
	char *person_missingness = "ex_0005.pm.dat"; // Person missingness data for unfiltered data-set
	char *founders = "ex_0005.ff.dat";	// Gender Data
	char *phenoName = "QT"; // Column to read from the phenotype file

	// Create the source data-set.
	// This is a one-off operation with the binary encoded data-set persisting somewhere on 
	// the file system.
	fill_binary_ped_command(command, map, ped, id, LDIR, RDIR);
	echo_command(command);
	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// First call of filter program to process and filter an allelic data-set.
	mode = 1;
	fill_pformat_data2_command(command, NUMPROCS, bmap, bped, phe, phenoName, 
							   NUMPROCS, run_id, LDIR, RDIR, rowMajor, nperm, mode, mind, geno, maf,
							   person_missingness, founders, hwe);

	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// 2nd call to join the fragments and create the filtered data-set.
	// This call is always a single core operation.
	mode = 2;
	fill_pformat_data2_command(command, 1, bmap, bped, phe, phenoName, 
							   NUMPROCS, run_id, LDIR, RDIR, rowMajor, nperm, mode, mind, geno, maf,
							   person_missingness, founders, hwe);

	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// On a 560K SNP data-set of 1000, 32 seconds wre required to loop through the 
	// data-set and a further 10 seconds to write the output file.
	// Running a filter job on 4 cores + improvements to the file I/O
	// dropped the wait time to 15 seconds for a filter job to complete.

	printf("\nFile Outputs:-\n");
	printf("\nrun_ex_0005.bmap.dat = Filtered SNP List\n");
	printf("run_ex_0005.fslog.txt = log file containing filter settings\n");
	printf("run_ex_0005.perm.dat = Permutation indices file\n");
	printf("run_ex_0005.ptd.dat = Persons to delete list for strata generation\n");
	printf("run_ex_0005.rseed.dat = Array of Random seeds for cores on a cluster\n");
	printf("run_ex_0005.rseed.txt = Initialisation value used by PRNG\n");
	printf("run_ex_0005.samples.txt = Sample array containg gender and phenotype data\n");
	printf("run_ex_0005.iped.dat = BitString array of allelic data (Individual major)\n");
	echo_spacer();
}

void ex_0006_association_scan_bt(void) {
	int result;
	char command[1000];

	// Text PLINK-compatible input files.
	char *map = "../data/snp_1000.map"; // Map file (SNP list)
	char *ped = "../data/snp_1000.ped"; // Pedigree File (Allelic Data)
	char *phe = "../data/random_phe_data.txt"; // Phenotype data

	// Filter paramaters
	char *rowMajor = "s"; // SNPs
	double mind = 0.03; // Missingness threshold person
	double geno = 0.031; // SNP missingness rate
	double maf = 0.15; // Major allelic frequency
	double hwe = 0.10; // HWE inclusion threshold
	int mode = 1; // Filtering Step (1 or 2)
	int nperm = 1000; // No. permutations required

	// Files of the binary encoded data-set created by "binary_ped.exe"
	char *id = "ex_0006"; // binary encoded data-set ID
	char *bmap = "ex_0006.bmap.dat"; // Data-set SNP list
	char *bped = "ex_0006.bped.dat"; // Data-set allelic data
	char *person_missingness = "ex_0006.pm.dat"; // Person missingness data for unfiltered data-set
	char *founders = "ex_0006.ff.dat";  // Gender Data
	char *phenoName = "BT";

	// Filtered data-set files
	char *run_id = "run_ex_0006"; // Filtered data-set ID
	char *samples = "run_ex_0006.samples.dat";
	char *filtered_bmap = "run_ex_0006.bmap.dat";
	char *filtered_sped = "run_ex_0006.sped.dat";
	char *permutation_indices = "run_ex_0006.perm.dat";

	// Create the source data-set.
	// This is a one-off operation with the binary encoded data-set persisting somewhere on 
	// the file system.
	fill_binary_ped_command(command, map, ped, id, LDIR, RDIR);
	echo_command(command);
	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// First call of filter program to process and filter an allelic data-set.
	mode = 1;
	fill_pformat_data2_command(command, NUMPROCS, bmap, bped, phe, phenoName, 
							   NUMPROCS, run_id, LDIR, RDIR, rowMajor, nperm, mode, mind, geno, maf,
							   person_missingness, founders, hwe);

	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// 2nd call to join the fragments and create the filtered data-set.
	// This call is always a single core operation.
	mode = 2;
	fill_pformat_data2_command(command, 1, bmap, bped, phe, phenoName, 
							   NUMPROCS, run_id, LDIR, RDIR, rowMajor, nperm, mode, mind, geno, maf,
							   person_missingness, founders, hwe);

	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// If reach here, a filtered data-set if ready to be shunted onto a
	// cluster for analysis and permutation.

	// The 1st call to 'passoc_bt_ft.exe' causes each processor nodes allocated to the calculation 
	// to run the non-permuted calculation and then uses the generated statistics to drive the permutation.
	mode = 1;
	fill_passoc_bt_ft_command(command, NUMPROCS, run_id, samples, filtered_bmap, filtered_sped, LDIR, RDIR, permutation_indices, mode);
	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// The 2nd call to 'passoc_bt_ft.exe' is to re-run the non-permuted calculation on a single core,
	// merge in the permutation results and create the final CSV data.
	// The 2nd call also uses the sort-by parameter to specify what column effects the sorting
	// of the final CSV file containing the statistical output.
	mode = 2;
	fill_passoc_bt_ft_command(command, 1, run_id, samples, filtered_bmap, filtered_sped, LDIR, RDIR, permutation_indices, mode);
	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	printf("\nFile Outputs:-\n");
	printf("run_ex_0006.assoc.csv = Sorted Results of the BT Association Scan\n");
	printf("run_ex_0006.fslog.txt = Log file containing filter settings\n");
	printf("run_ex_0006.rseed.txt = Initialisation value used by PRNG\n");
	echo_spacer(); 
}


void ex_0007_association_scan_bt_ft(void) {
	int result;
	char command[1000];

	// Text PLINK-compatible input files.
	char *map = "../data/snp_1000.map"; // Map file (SNP list)
	char *ped = "../data/snp_1000.ped"; // Pedigree File (Allelic Data)
	char *phe = "../data/random_phe_data.txt"; // Phenotype data

	// Filter paramaters
	char *rowMajor = "s"; // SNPs
	double mind = 0.03; // Missingness threshold person
	double geno = 0.031; // SNP missingness rate
	double maf = 0.15; // Major allelic frequency
	double hwe = 0.10; // HWE inclusion threshold
	int mode = 1; // Filtering Step (1 or 2)
	int nperm = 1000; // No. permutations required

	// Files of the binary encoded data-set created by "binary_ped.exe"
	char *id = "ex_0007"; // binary encoded data-set ID
	char *bmap = "ex_0007.bmap.dat"; // Data-set SNP list
	char *bped = "ex_0007.bped.dat"; // Data-set allelic data
	char *person_missingness = "ex_0007.pm.dat"; // Person missingness data for unfiltered data-set
	char *founders = "ex_0007.ff.dat";  // Gender Data
	char *phenoName = "BT";

	// Filtered data-set files
	char *run_id = "run_ex_0007"; // Filtered data-set ID
	char *samples = "run_ex_0007.samples.dat";
	char *filtered_bmap = "run_ex_0007.bmap.dat";
	char *filtered_sped = "run_ex_0007.sped.dat";
	char *permutation_indices = "run_ex_0007.perm.dat";

	// Create the source data-set.
	// This is a one-off operation with the binary encoded data-set persisting somewhere on 
	// the file system.
	fill_binary_ped_command(command, map, ped, id, LDIR, RDIR);
	echo_command(command);
	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// First call of filter program to process and filter an allelic data-set.
	mode = 1;
	fill_pformat_data2_command(command, NUMPROCS, bmap, bped, phe, phenoName, 
							   NUMPROCS, run_id, LDIR, RDIR, rowMajor, nperm, mode, mind, geno, maf,
							   person_missingness, founders, hwe);

	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// 2nd call to join the fragments and create the filtered data-set.
	// This call is always a single core operation.
	mode = 2;
	fill_pformat_data2_command(command, 1, bmap, bped, phe, phenoName, 
							   NUMPROCS, run_id, LDIR, RDIR, rowMajor, nperm, mode, mind, geno, maf,
							   person_missingness, founders, hwe);

	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// If reach here, a filtered data-set if ready to be shunted onto a
	// cluster for analysis and permutation.

	// The 1st call to 'passoc_bt_ft.exe' causes each processor nodes allocated to the calculation 
	// to run the non-permuted calculation and then uses the generated statistics to drive the permutation.
	mode = 1;
	fill_passoc_bt_ft_command(command, NUMPROCS, run_id, samples, filtered_bmap, filtered_sped, LDIR, RDIR, permutation_indices, mode);
	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// This call inserts parameters to the 'passoc_bt_ft.exe' so that a single core
	// can repeat a section of a calculation which failed on a node (computers do break occassionally).
	// The node failure is on rank=0 of 2 processors allocated to the calculation.
	mode = 3;
	fill_passoc_bt_ft_command(command, 1, run_id, samples, filtered_bmap, filtered_sped, LDIR, RDIR, permutation_indices, mode);
	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;
	
	// The 3rd call to 'passoc_bt_ft.exe' is to re-run the non-permuted calculation on a single core,
	// merge in the permutation results and create the final CSV data.
	// The 3rd call also uses the sort-by parameter to specify what column effects the sorting
	// of the final CSV file containing the statistical output.
	mode = 2;
	fill_passoc_bt_ft_command(command, 1, run_id, samples, filtered_bmap, filtered_sped, LDIR, RDIR, permutation_indices, mode);
	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	printf("\nFile Outputs:-\n");
	printf("run_ex_0007.assoc.csv = Sorted Results of the BT Association Scan\n");
	printf("run_ex_0007.fslog.txt = Log file containing filter settings\n");
	printf("run_ex_0007.rseed.txt = Initialisation value used by PRNG\n");
	echo_spacer(); 
}


void ex_0008_association_scan_qt_ft(void) {
	int result;
	char command[1000];

	// Text PLINK-compatible input files.
	char *map = "../data/snp_1000.map"; // Map file (SNP list)
	char *ped = "../data/snp_1000.ped"; // Pedigree File (Allelic Data)
	char *phe = "../data/random_phe_data.txt"; // Phenotype data

	// Filter paramaters
	char *rowMajor = "s"; // SNPs
	double mind = 0.03; // Missingness threshold person
	double geno = 0.031; // SNP missingness rate
	double maf = 0.15; // Major allelic frequency
	double hwe = 0.10; // HWE inclusion threshold
	int mode = 1; // Filtering Step (1 or 2)
	int nperm = 1000; // No. permutations required

	// Files of the binary encoded data-set created by "binary_ped.exe"
	char *id = "ex_0008"; // binary encoded data-set ID
	char *bmap = "ex_0008.bmap.dat"; // Data-set SNP list
	char *bped = "ex_0008.bped.dat"; // Data-set allelic data
	char *person_missingness = "ex_0008.pm.dat"; // Person missingness data for unfiltered data-set
	char *founders = "ex_0008.ff.dat";  // Gender Data
	char *phenoName = "QT";

	// Filtered data-set files
	char *run_id = "run_ex_0008"; // Filtered data-set ID
	char *samples = "run_ex_0008.samples.dat";
	char *filtered_bmap = "run_ex_0008.bmap.dat";
	char *filtered_sped = "run_ex_0008.sped.dat";
	char *permutation_indices = "run_ex_0008.perm.dat";

	// Create the source data-set.
	// This is a one-off operation with the binary encoded data-set persisting somewhere on 
	// the file system.
	fill_binary_ped_command(command, map, ped, id, LDIR, RDIR);
	echo_command(command);
	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// First call of filter program to process and filter an allelic data-set.
	mode = 1;
	fill_pformat_data2_command(command, NUMPROCS, bmap, bped, phe, phenoName, 
							   NUMPROCS, run_id, LDIR, RDIR, rowMajor, nperm, mode, mind, geno, maf,
							   person_missingness, founders, hwe);

	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// 2nd call to join the fragments and create the filtered data-set.
	// This call is always a single core operation.
	mode = 2;
	fill_pformat_data2_command(command, 1, bmap, bped, phe, phenoName, 
							   NUMPROCS, run_id, LDIR, RDIR, rowMajor, nperm, mode, mind, geno, maf,
							   person_missingness, founders, hwe);

	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// If reach here, a filtered data-set if ready to be shunted onto a
	// cluster for analysis and permutation.

	// The 1st call to 'passoc_qt_ft.exe' causes each processor nodes allocated to the calculation 
	// to run the non-permuted calculation and then uses the generated statistics to drive the permutation.
	mode = 1;
	fill_passoc_qt_ft_command(command, NUMPROCS, run_id, samples, filtered_bmap, filtered_sped, LDIR, RDIR, permutation_indices, mode);
	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// This call inserts parameters to the 'passoc_qt_ft.exe' so that a single core
	// can repeat a section of a calculation which failed on a node (computers do break occassionally).
	// The node failure is on rank=0 of 2 processors allocated to the calculation.
	mode = 3;
	fill_passoc_qt_ft_command(command, 1, run_id, samples, filtered_bmap, filtered_sped, LDIR, RDIR, permutation_indices, mode);
	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;
	
	// The 3rd call to 'passoc_bt_ft.exe' is to re-run the non-permuted calculation on a single core,
	// merge in the permutation results and create the final CSV data.
	// The 3rd call also uses the sort-by parameter to specify what column effects the sorting
	// of the final CSV file containing the statistical output.
	mode = 2;
	fill_passoc_qt_ft_command(command, 1, run_id, samples, filtered_bmap, filtered_sped, LDIR, RDIR, permutation_indices, mode);
	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	printf("\nFile Outputs:-\n");
	printf("run_ex_0008.assoc.csv = Sorted Results of the QT Association Scan\n");
	printf("run_ex_0008.fslog.txt = Log file containing filter settings\n");
	printf("run_ex_0008.rseed.txt = Initialisation value used by PRNG\n");
	echo_spacer(); 
}


void ex_0009_model_association_ft(void) {
	int result;
	char command[1000];

	// Text PLINK-compatible input files.
	char *map = "../data/snp_1000.map"; // Map file (SNP list)
	char *ped = "../data/snp_1000.ped"; // Pedigree File (Allelic Data)
	char *phe = "../data/random_phe_data.txt"; // Phenotype data

	// Filter paramaters
	char *rowMajor = "s"; // SNPs
	double mind = 0.03; // Missingness threshold person
	double geno = 0.031; // SNP missingness rate
	double maf = 0.15; // Major allelic frequency
	double hwe = 0.10; // HWE inclusion threshold
	int mode = 1; // Filtering Step (1 or 2)
	int nperm = 1000; // No. permutations required

	// Files of the binary encoded data-set created by "binary_ped.exe"
	char *id = "ex_0009"; // binary encoded data-set ID
	char *bmap = "ex_0009.bmap.dat"; // Data-set SNP list
	char *bped = "ex_0009.bped.dat"; // Data-set allelic data
	char *person_missingness = "ex_0009.pm.dat"; // Person missingness data for unfiltered data-set
	char *founders = "ex_0009.ff.dat";  // Gender Data
	char *phenoName = "BT";

	// Filtered data-set files
	char *run_id = "run_ex_0009"; // Filtered data-set ID
	char *samples = "run_ex_0009.samples.dat";
	char *filtered_bmap = "run_ex_0009.bmap.dat";
	char *filtered_sped = "run_ex_0009.sped.dat";
	char *permutation_indices = "run_ex_0009.perm.dat";

	// Create the source data-set.
	// This is a one-off operation with the binary encoded data-set persisting somewhere on 
	// the file system.
	fill_binary_ped_command(command, map, ped, id, LDIR, RDIR);
	echo_command(command);
	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// First call of filter program to process and filter an allelic data-set.
	mode = 1;
	fill_pformat_data2_command(command, NUMPROCS, bmap, bped, phe, phenoName, 
							   NUMPROCS, run_id, LDIR, RDIR, rowMajor, nperm, mode, mind, geno, maf,
							   person_missingness, founders, hwe);

	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// 2nd call to join the fragments and create the filtered data-set.
	// This call is always a single core operation.
	mode = 2;
	fill_pformat_data2_command(command, 1, bmap, bped, phe, phenoName, 
							   NUMPROCS, run_id, LDIR, RDIR, rowMajor, nperm, mode, mind, geno, maf,
							   person_missingness, founders, hwe);

	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// If reach here, a filtered data-set if ready to be shunted onto a
	// cluster for analysis and permutation.

	// The 1st call to 'pmodel_ft.exe' causes each processor nodes allocated to the calculation 
	// to run the non-permuted calculation and then uses the generated statistics to drive the permutation.
	mode = 1;
	fill_pmodel_ft_command(command, NUMPROCS, run_id, samples, filtered_bmap, filtered_sped, LDIR, RDIR, permutation_indices, mode);
	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// This call inserts parameters to the 'pmodel_ft.exe' so that a single core
	// can repeat a section of a calculation which failed on a node (computers do break occassionally).
	// The node failure is on rank=0 of 2 processors allocated to the calculation.
	mode = 3;
	fill_pmodel_ft_command(command, 1, run_id, samples, filtered_bmap, filtered_sped, LDIR, RDIR, permutation_indices, mode);
	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;
	
	// The 3rd call to 'passoc_bt_ft.exe' is to re-run the non-permuted calculation on a single core,
	// merge in the permutation results and create the final CSV data.
	// The 3rd call also uses the sort-by parameter to specify what column effects the sorting
	// of the final CSV file containing the statistical output.
	mode = 2;
	fill_pmodel_ft_command(command, 1, run_id, samples, filtered_bmap, filtered_sped, LDIR, RDIR, permutation_indices, mode);
	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	printf("\nFile Outputs:-\n");
	printf("run_ex_0009.model.csv = Sorted Results of the Model Association Scan\n");
	printf("run_ex_0009.fslog.txt = Log file containing filter settings\n");
	printf("run_ex_0009.rseed.txt = Initialisation value used by PRNG\n");
	echo_spacer();
}

void ex_0010_fast_epistasis_bt(void) {
	int result;
	char command[1000];

	// Text PLINK-compatible input files.
	char *map = "../data/snp_100.map"; // Map file (SNP list)
	char *ped = "../data/snp_100.ped"; // Pedigree File (Allelic Data)
	char *phe = "../data/random_phe_data.txt"; // Phenotype data

	// Filter paramaters
	char *rowMajor = "s"; // SNPs
	double mind = 0.03; // Missingness threshold person
	double geno = 0.031; // SNP missingness rate
	double maf = 0.15; // Major allelic frequency
	double hwe = 0.10; // HWE inclusion threshold
	int mode = 1; // Filtering Step (1 or 2)
	int nperm = 2; // No. permutations required, low as interaction calculations heavy operation

	// Files of the binary encoded data-set created by "binary_ped.exe"
	char *id = "ex_0010"; // binary encoded data-set ID
	char *bmap = "ex_0010.bmap.dat"; // Data-set SNP list
	char *bped = "ex_0010.bped.dat"; // Data-set allelic data
	char *person_missingness = "ex_0010.pm.dat"; // Person missingness data for unfiltered data-set
	char *founders = "ex_0010.ff.dat";  // Gender Data
	char *phenoName = "BT";

	// Filtered data-set files
	char *run_id = "run_ex_0010"; // Filtered data-set ID
	char *samples = "run_ex_0010.samples.dat";
	char *filtered_bmap = "run_ex_0010.bmap.dat";
	char *filtered_sped = "run_ex_0010.sped.dat";

	// Create the source data-set.
	// This is a one-off operation with the binary encoded data-set persisting somewhere on 
	// the file system.
	fill_binary_ped_command(command, map, ped, id, LDIR, RDIR);
	echo_command(command);
	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// First call of filter program to process and filter an allelic data-set.
	mode = 1;
	fill_pformat_data2_command(command, NUMPROCS, bmap, bped, phe, phenoName, 
							   NUMPROCS, run_id, LDIR, RDIR, rowMajor, nperm, mode, mind, geno, maf,
							   person_missingness, founders, hwe);

	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// 2nd call to join the fragments and create the filtered data-set.
	// This call is always a single core operation.
	mode = 2;
	fill_pformat_data2_command(command, 1, bmap, bped, phe, phenoName, 
							   NUMPROCS, run_id, LDIR, RDIR, rowMajor, nperm, mode, mind, geno, maf,
							   person_missingness, founders, hwe);

	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// First call to 'fast_pepi_bt.exe', no permuted calculation spread across N-cores.
	mode = MODE_NON_PERMUTED;
	fill_fast_pepi_bt_v2_command(command, NUMPROCS, run_id, samples, filtered_bmap, filtered_sped, LDIR, RDIR, nperm, mode);
	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// 2nd call on a single core is to combine the results from all ranks and apply a p-value threshold to find
	// putative interesting interactions in the data-set.
	mode = MODE_NON_PERMUTED_MERGE;
	fill_fast_pepi_bt_v2_command(command, 1, run_id, samples, filtered_bmap, filtered_sped, LDIR, RDIR, nperm, mode);
	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// 3rd call on a single core to write out a final output, sorted by P-value.
	// The 2nd call is shared by this non-permuted calculation and the permuted calcualtion.
	// This functional sharing explains this 3-stage invocation for generating non-permuted results.
	mode = MODE_STOP;
	fill_fast_pepi_bt_v2_command(command, 1, run_id, samples, filtered_bmap, filtered_sped, LDIR, RDIR, nperm, mode);
	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// As the code is not doing permutation, the randomised result data columns (EMP1 and EMP2)
	// are set as one.
	printf("\nFile Outputs:-\n");
	printf("\nrun_ex_0010.fepi.csv = Sorted list of Interesting Interactions\n");
	printf("run_ex_0010.fslog.txt = log file containing filter settings\n");
	printf("run_ex_0010.rseed.txt = Initialisation value used by PRNG\n");
	echo_spacer(); 
}

void ex_0011_fast_epistasis_bt_with_permutation(void) {
	int result;
	char command[1000];

	// Text PLINK-compatible input files.
	char *map = "../data/snp_100.map"; // Map file (SNP list)
	char *ped = "../data/snp_100.ped"; // Pedigree File (Allelic Data)
	char *phe = "../data/random_phe_data.txt"; // Phenotype data

	// Filter paramaters
	char *rowMajor = "s"; // SNPs
	double mind = 0.03; // Missingness threshold person
	double geno = 0.031; // SNP missingness rate
	double maf = 0.15; // Major allelic frequency
	double hwe = 0.10; // HWE inclusion threshold
	int mode = 1; // Filtering Step (1 or 2)
	int nperm = 100; // No. permutations required, low as interaction calculations heavy operation

	// Files of the binary encoded data-set created by "binary_ped.exe"
	char *id = "ex_0011"; // binary encoded data-set ID
	char *bmap = "ex_0011.bmap.dat"; // Data-set SNP list
	char *bped = "ex_0011.bped.dat"; // Data-set allelic data
	char *person_missingness = "ex_0011.pm.dat"; // Person missingness data for unfiltered data-set
	char *founders = "ex_0011.ff.dat";  // Gender Data
	char *phenoName = "BT";

	// Filtered data-set files
	char *run_id = "run_ex_0011"; // Filtered data-set ID
	char *samples = "run_ex_0011.samples.dat";
	char *filtered_bmap = "run_ex_0011.bmap.dat";
	char *filtered_sped = "run_ex_0011.sped.dat";

	// Create the source data-set.
	// This is a one-off operation with the binary encoded data-set persisting somewhere on 
	// the file system.
	fill_binary_ped_command(command, map, ped, id, LDIR, RDIR);
	echo_command(command);
	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// First call of filter program to process and filter an allelic data-set.
	mode = 1;
	fill_pformat_data2_command(command, NUMPROCS, bmap, bped, phe, phenoName, 
							   NUMPROCS, run_id, LDIR, RDIR, rowMajor, nperm, mode, mind, geno, maf,
							   person_missingness, founders, hwe);

	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// 2nd call to join the fragments and create the filtered data-set.
	// This call is always a single core operation.
	mode = 2;
	fill_pformat_data2_command(command, 1, bmap, bped, phe, phenoName, 
							   NUMPROCS, run_id, LDIR, RDIR, rowMajor, nperm, mode, mind, geno, maf,
							   person_missingness, founders, hwe);

	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// First call to 'fast_pepi_bt.exe', no permuted calculation spread across N-cores.
	mode = MODE_NON_PERMUTED;
	fill_fast_pepi_bt_v2_command(command, NUMPROCS, run_id, samples, filtered_bmap, filtered_sped, LDIR, RDIR, nperm, mode);
	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// 2nd call on a single core is to combine the results from all ranks and apply a p-value threshold to find
	// putative interesting interactions in the data-set.
	mode = MODE_NON_PERMUTED_MERGE;
	fill_fast_pepi_bt_v2_command(command, 1, run_id, samples, filtered_bmap, filtered_sped, LDIR, RDIR, nperm, mode);
	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// 3rd call on N-cores to run the permutation.
	// each core runs a randomised data-set so running epistasis on large data-set is
	// not advised as it will take ages.
	mode = MODE_PERMUTED;
	fill_fast_pepi_bt_v2_command(command, NUMPROCS, run_id, samples, filtered_bmap, filtered_sped, LDIR, RDIR, nperm, mode);
	printf("\n");
	result = system(command);
	printf("\n");
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// 4th call on single core to merge the permutation results into a single file.
	// Single core operation
	mode = MODE_PERMUTED_MERGE;
	fill_fast_pepi_bt_v2_command(command, 1, run_id, samples, filtered_bmap, filtered_sped, LDIR, RDIR, nperm, mode);
	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// As the code is not doing permutation, the randomised result data columns (EMP1 and EMP2)
	// are set as zero.
	printf("\nFile Outputs:-\n");
	printf("\nrun_ex_0011.fepi.csv = Sorted list of Interesting Interactions\n");
	printf("run_ex_0011.fslog.txt = log file containing filter settings\n");
	printf("run_ex_0011.rseed.txt = Initialisation value used by PRNG\n");
	echo_spacer(); 
}

void ex_0012_stratified_assoc_scan_bt(void) {
	int result;
	char command[1000];

	// Text PLINK-compatible input files.
	char *map = "../data/snp_1000.map"; // Map file (SNP list)
	char *ped = "../data/snp_1000.ped"; // Pedigree File (Allelic Data)
	char *phe = "../data/random_phe_data.txt"; // Phenotype data

	// Filter paramaters
	char *rowMajor = "s"; // SNPs
	double mind = 0.03; // Missingness threshold person
	double geno = 0.031; // SNP missingness rate
	double maf = 0.15; // Major allelic frequency
	double hwe = 0.10; // HWE inclusion threshold
	int mode = 1; // Filtering Step (1 or 2)
	int nperm = 2; // No. permutations required

	// Files of the binary encoded data-set created by "binary_ped.exe"
	char *id = "ex_0012"; // binary encoded data-set ID
	char *bmap = "ex_0012.bmap.dat"; // Data-set SNP list
	char *bped = "ex_0012.bped.dat"; // Data-set allelic data
	char *person_missingness = "ex_0012.pm.dat"; // Person missingness data for unfiltered data-set
	char *founders = "ex_0012.ff.dat";  // Gender Data
	char *phenoName = "BT";
	char *strataPhenoName = "COUNTRY1"; // Phenotype column to generate the strata information.

	// Filtered data-set files
	char *run_id = "run_ex_0012"; // Filtered data-set ID
	char *samples = "run_ex_0012.samples.dat";
	char *filtered_bmap = "run_ex_0012.bmap.dat";
	char *filtered_sped = "run_ex_0012.sped.dat";
	char *permutation_indices = "run_ex_0012.perm.dat"; // Stratified permutation index file
	char *strataFile = "run_ex_0012.sol.dat";
	double ci = 0.95; // Confidence Interval

	// Create the source data-set.
	// This is a one-off operation with the binary encoded data-set persisting somewhere on 
	// the file system.
	fill_binary_ped_command(command, map, ped, id, LDIR, RDIR);
	echo_command(command);
	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// First call of filter program to process and filter an allelic data-set.
	mode = 1;
	fill_pformat_data2_command(command, NUMPROCS, bmap, bped, phe, phenoName, 
							   NUMPROCS, run_id, LDIR, RDIR, rowMajor, nperm, mode, mind, geno, maf,
							   person_missingness, founders, hwe);

	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// 2nd call to join the fragments and create the filtered data-set.
	// This call is always a single core operation.
	mode = 2;
	fill_pformat_data2_command(command, 1, bmap, bped, phe, phenoName, 
							   NUMPROCS, run_id, LDIR, RDIR, rowMajor, nperm, mode, mind, geno, maf,
							   person_missingness, founders, hwe);

	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// Call 'solcol.exe' on a single core to create a strata membership file.
	// Solcol must be called after pformat in order to generate data the PTD file
	// listing what people were removed from a data-set due to filtering.
	// 'Solcol.exe' generates a stratified permutation file to control
	// stratified randomisation on a cluster, replacing output from 'pformat_data2.exe'.
	// To suppress randomisation artifacts being created by 'pformat_data2.exe', use the '--no-rng'
	// at the command-line.
	fill_solcol_command(command, 1, run_id, phe, strataPhenoName, LDIR, RDIR, nperm);
	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// If reach here, a filtered data-set if ready to be shunted onto a
	// cluster for analysis and permutation.

	// The 1st call to 'pcmh_bt.exe' causes each processor nodes allocated to the calculation 
	// to run the non-permuted calculation and then uses the generated statistics to drive the permutation.
	mode = 1;
	fill_pcmh_ft_command(command, NUMPROCS, run_id, samples, filtered_bmap, filtered_sped, LDIR, RDIR, permutation_indices, mode, strataFile, ci);
	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// The 2nd call to 'passoc_bt_ft.exe' is to re-run the non-permuted calculation on a single core,
	// merge in the permutation results and create the final CSV data.
	// The 2nd call also uses the sort-by parameter to specify what column effects the sorting
	// of the final CSV file containing the statistical output.
	mode = 2;
	fill_pcmh_ft_command(command, 1, run_id, samples, filtered_bmap, filtered_sped, LDIR, RDIR, permutation_indices, mode, strataFile, ci);
	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	printf("\nFile Outputs:-\n");
	printf("run_ex_0012.cmh.csv = Sorted Results of the CMH Stratified Association Scan\n");
	printf("run_ex_0012.fslog.txt = Log file containing filter settings\n");
	printf("run_ex_0012.rseed.txt = Initialisation value used by PRNG\n");
	echo_spacer();
}

void ex_0013_slow_epistasis_bt(void) {
	int result;
	char command[1000];

	// Text PLINK-compatible input files.
	char *map = "../data/snp_100.map"; // Map file (SNP list)
	char *ped = "../data/snp_100.ped"; // Pedigree File (Allelic Data)
	char *phe = "../data/random_phe_data.txt"; // Phenotype data

	// Filter paramaters
	char *rowMajor = "i"; // SNPs
	double mind = 0.03; // Missingness threshold person
	double geno = 0.031; // SNP missingness rate
	double maf = 0.15; // Major allelic frequency
	double hwe = 0.10; // HWE inclusion threshold
	int mode = 1; // Filtering Step (1 or 2)
	int nperm = 2; // No. permutations required, low as interaction calculations heavy operation

	// Files of the binary encoded data-set created by "binary_ped.exe"
	char *id = "ex_0013"; // binary encoded data-set ID
	char *bmap = "ex_0013.bmap.dat"; // Data-set SNP list
	char *bped = "ex_0013.bped.dat"; // Data-set allelic data
	char *person_missingness = "ex_0013.pm.dat"; // Person missingness data for unfiltered data-set
	char *founders = "ex_0013.ff.dat";  // Gender Data
	char *phenoName = "BT";

	// Filtered data-set files
	char *run_id = "run_ex_0013"; // Filtered data-set ID
	char *samples = "run_ex_0013.samples.dat";
	char *filtered_bmap = "run_ex_0013.bmap.dat";
	char *filtered_iped = "run_ex_0013.iped.dat";

	// Create the source data-set.
	// This is a one-off operation with the binary encoded data-set persisting somewhere on 
	// the file system.
	fill_binary_ped_command(command, map, ped, id, LDIR, RDIR);
	echo_command(command);
	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// First call of filter program to process and filter an allelic data-set.
	mode = 1;
	fill_pformat_data2_command(command, NUMPROCS, bmap, bped, phe, phenoName, 
							   NUMPROCS, run_id, LDIR, RDIR, rowMajor, nperm, mode, mind, geno, maf,
							   person_missingness, founders, hwe);

	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// 2nd call to join the fragments and create the filtered data-set.
	// This call is always a single core operation.
	mode = 2;
	fill_pformat_data2_command(command, 1, bmap, bped, phe, phenoName, 
							   NUMPROCS, run_id, LDIR, RDIR, rowMajor, nperm, mode, mind, geno, maf,
							   person_missingness, founders, hwe);

	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// First call to 'pepi_bt_v2.exe', no permuted calculation spread across N-cores.
	mode = MODE_NON_PERMUTED;
	fill_pepi_bt_v2_command(command, NUMPROCS, run_id, samples, filtered_bmap, filtered_iped, LDIR, RDIR, nperm, mode);
	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// 2nd call on a single core is to combine the results from all ranks and apply a p-value threshold to find
	// putative interesting interactions in the data-set.
	mode = MODE_NON_PERMUTED_MERGE;
	fill_pepi_bt_v2_command(command, 1, run_id, samples, filtered_bmap, filtered_iped, LDIR, RDIR, nperm, mode);
	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// 3rd call on a single core to write out a final output, sorted by P-value.
	// The 2nd call is shared by this non-permuted calculation and the permuted calcualtion.
	// This functional sharing explains this 3-stage invocation for generating non-permuted results.
	mode = MODE_STOP;
	fill_pepi_bt_v2_command(command, 1, run_id, samples, filtered_bmap, filtered_iped, LDIR, RDIR, nperm, mode);
	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	printf("\nFile Outputs:-\n");
	printf("\nrun_ex_0013.fepi.csv = Interesting Interactions for BT Slow Epistasis\n");
	printf("run_ex_0013.fslog.txt = log file containing filter settings\n");
	printf("run_ex_0013.rseed.txt = Initialisation value used by PRNG\n");
	echo_spacer(); 
}


void ex_0014_slow_epistasis_bt_with_permutation(void) {
	int result;
	char command[1000];

	// Text PLINK-compatible input files.
	char *map = "../data/snp_100.map"; // Map file (SNP list)
	char *ped = "../data/snp_100.ped"; // Pedigree File (Allelic Data)
	char *phe = "../data/random_phe_data.txt"; // Phenotype data

	// Filter paramaters
	char *rowMajor = "i"; // SNPs
	double mind = 0.03; // Missingness threshold person
	double geno = 0.031; // SNP missingness rate
	double maf = 0.15; // Major allelic frequency
	double hwe = 0.10; // HWE inclusion threshold
	int mode = 1; // Filtering Step (1 or 2)
	int nperm = 100; // No. permutations required, low as interaction calculations heavy operation

	// Files of the binary encoded data-set created by "binary_ped.exe"
	char *id = "ex_0014"; // binary encoded data-set ID
	char *bmap = "ex_0014.bmap.dat"; // Data-set SNP list
	char *bped = "ex_0014.bped.dat"; // Data-set allelic data
	char *person_missingness = "ex_0014.pm.dat"; // Person missingness data for unfiltered data-set
	char *founders = "ex_0014.ff.dat";  // Gender Data
	char *phenoName = "BT";

	// Filtered data-set files
	char *run_id = "run_ex_0014"; // Filtered data-set ID
	char *samples = "run_ex_0014.samples.dat";
	char *filtered_bmap = "run_ex_0014.bmap.dat";
	char *filtered_sped = "run_ex_0014.iped.dat";

	// Create the source data-set.
	// This is a one-off operation with the binary encoded data-set persisting somewhere on 
	// the file system.
	fill_binary_ped_command(command, map, ped, id, LDIR, RDIR);
	echo_command(command);
	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// First call of filter program to process and filter an allelic data-set.
	mode = 1;
	fill_pformat_data2_command(command, NUMPROCS, bmap, bped, phe, phenoName, 
							   NUMPROCS, run_id, LDIR, RDIR, rowMajor, nperm, mode, mind, geno, maf,
							   person_missingness, founders, hwe);

	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// 2nd call to join the fragments and create the filtered data-set.
	// This call is always a single core operation.
	mode = 2;
	fill_pformat_data2_command(command, 1, bmap, bped, phe, phenoName, 
							   NUMPROCS, run_id, LDIR, RDIR, rowMajor, nperm, mode, mind, geno, maf,
							   person_missingness, founders, hwe);

	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// First call to 'pepi_bt_v2.exe', no permuted calculation spread across N-cores.
	mode = MODE_NON_PERMUTED;
	fill_pepi_bt_v2_command(command, NUMPROCS, run_id, samples, filtered_bmap, filtered_sped, LDIR, RDIR, nperm, mode);
	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// 2nd call on a single core is to combine the results from all ranks and apply a p-value threshold to find
	// putative interesting interactions in the data-set.
	mode = MODE_NON_PERMUTED_MERGE;
	fill_pepi_bt_v2_command(command, 1, run_id, samples, filtered_bmap, filtered_sped, LDIR, RDIR, nperm, mode);
	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// 3rd call on N-cores to run the permutation.
	// each core runs a randomised data-set so running epistasis on large data-set is
	// not advised as it will take ages.
	mode = MODE_PERMUTED;
	fill_pepi_bt_v2_command(command, NUMPROCS, run_id, samples, filtered_bmap, filtered_sped, LDIR, RDIR, nperm, mode);
	printf("\n");
	result = system(command);
	printf("\n");
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// 4th call on single core to merge the permutation results into a single file.
	// Single core operation
	mode = MODE_PERMUTED_MERGE;
	fill_pepi_bt_v2_command(command, 1, run_id, samples, filtered_bmap, filtered_sped, LDIR, RDIR, nperm, mode);
	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// As the code is not doing permutation, the randomised result data columns (EMP1 and EMP2)
	// are set as zero.
	printf("\nFile Outputs:-\n");
	printf("\nrun_ex_0014.fepi.csv = Sorted list of Interesting Interactions with permutation statistics\n");
	printf("run_ex_0014.fslog.txt = log file containing filter settings\n");
	printf("run_ex_0014.rseed.txt = Initialisation value used by PRNG\n");
	echo_spacer(); 
}

void ex_0015_slow_epistasis_qt(void) {
	int result;
	char command[1000];

	// Text PLINK-compatible input files.
	char *map = "../data/snp_100.map"; // Map file (SNP list)
	char *ped = "../data/snp_100.ped"; // Pedigree File (Allelic Data)
	char *phe = "../data/random_phe_data.txt"; // Phenotype data

	// Filter paramaters
	char *rowMajor = "i"; // SNPs
	double mind = 0.03; // Missingness threshold person
	double geno = 0.031; // SNP missingness rate
	double maf = 0.15; // Major allelic frequency
	double hwe = 0.10; // HWE inclusion threshold
	int mode = 1; // Filtering Step (1 or 2)
	int nperm = 2; // No. permutations required, low as interaction calculations heavy operation

	// Files of the binary encoded data-set created by "binary_ped.exe"
	char *id = "ex_0015"; // binary encoded data-set ID
	char *bmap = "ex_0015.bmap.dat"; // Data-set SNP list
	char *bped = "ex_0015.bped.dat"; // Data-set allelic data
	char *person_missingness = "ex_0015.pm.dat"; // Person missingness data for unfiltered data-set
	char *founders = "ex_0015.ff.dat";  // Gender Data
	char *phenoName = "QT";

	// Filtered data-set files
	char *run_id = "run_ex_0015"; // Filtered data-set ID
	char *samples = "run_ex_0015.samples.dat";
	char *filtered_bmap = "run_ex_0015.bmap.dat";
	char *filtered_iped = "run_ex_0015.iped.dat";

	// Create the source data-set.
	// This is a one-off operation with the binary encoded data-set persisting somewhere on 
	// the file system.
	fill_binary_ped_command(command, map, ped, id, LDIR, RDIR);
	echo_command(command);
	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// First call of filter program to process and filter an allelic data-set.
	mode = 1;
	fill_pformat_data2_command(command, NUMPROCS, bmap, bped, phe, phenoName, 
							   NUMPROCS, run_id, LDIR, RDIR, rowMajor, nperm, mode, mind, geno, maf,
							   person_missingness, founders, hwe);

	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// 2nd call to join the fragments and create the filtered data-set.
	// This call is always a single core operation.
	mode = 2;
	fill_pformat_data2_command(command, 1, bmap, bped, phe, phenoName, 
							   NUMPROCS, run_id, LDIR, RDIR, rowMajor, nperm, mode, mind, geno, maf,
							   person_missingness, founders, hwe);

	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// First call to 'pepi_qt_v2.exe', non-permuted calculation spread across N-cores.
	mode = MODE_NON_PERMUTED;
	fill_pepi_qt_v2_command(command, NUMPROCS, run_id, samples, filtered_bmap, filtered_iped, LDIR, RDIR, nperm, mode);
	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// 2nd call on a single core is to combine the results from all ranks and apply a p-value threshold to find
	// putative interesting interactions in the data-set.
	mode = MODE_NON_PERMUTED_MERGE;
	fill_pepi_qt_v2_command(command, 1, run_id, samples, filtered_bmap, filtered_iped, LDIR, RDIR, nperm, mode);
	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// 3rd call on a single core to write out a final output, sorted by P-value.
	// The 2nd call is shared by this non-permuted calculation and the permuted calcualtion.
	// This functional sharing explains this 3-stage invocation for generating non-permuted results.
	mode = MODE_STOP;
	fill_pepi_qt_v2_command(command, 1, run_id, samples, filtered_bmap, filtered_iped, LDIR, RDIR, nperm, mode);
	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	printf("\nFile Outputs:-\n");
	printf("\nrun_ex_0015.fepi.csv = Interesting Interactions for QT Slow Epistasis\n");
	printf("run_ex_0015.fslog.txt = log file containing filter settings\n");
	printf("run_ex_0015.rseed.txt = Initialisation value used by PRNG\n");
	echo_spacer(); 
}

void ex_0016_slow_epistasis_qt_with_permutation(void) {
	int result;
	char command[1000];

	// Text PLINK-compatible input files.
	char *map = "../data/snp_100.map"; // Map file (SNP list)
	char *ped = "../data/snp_100.ped"; // Pedigree File (Allelic Data)
	char *phe = "../data/random_phe_data.txt"; // Phenotype data

	// Filter paramaters
	char *rowMajor = "i"; // SNPs
	double mind = 0.03; // Missingness threshold person
	double geno = 0.031; // SNP missingness rate
	double maf = 0.15; // Major allelic frequency
	double hwe = 0.10; // HWE inclusion threshold
	int mode = 1; // Filtering Step (1 or 2)
	int nperm = 100; // No. permutations required, low as interaction calculations heavy operation

	// Files of the binary encoded data-set created by "binary_ped.exe"
	char *id = "ex_0016"; // binary encoded data-set ID
	char *bmap = "ex_0016.bmap.dat"; // Data-set SNP list
	char *bped = "ex_0016.bped.dat"; // Data-set allelic data
	char *person_missingness = "ex_0016.pm.dat"; // Person missingness data for unfiltered data-set
	char *founders = "ex_0016.ff.dat";  // Gender Data
	char *phenoName = "QT";

	// Filtered data-set files
	char *run_id = "run_ex_0016"; // Filtered data-set ID
	char *samples = "run_ex_0016.samples.dat";
	char *filtered_bmap = "run_ex_0016.bmap.dat";
	char *filtered_sped = "run_ex_0016.iped.dat";

	// Create the source data-set.
	// This is a one-off operation with the binary encoded data-set persisting somewhere on 
	// the file system.
	fill_binary_ped_command(command, map, ped, id, LDIR, RDIR);
	echo_command(command);
	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// First call of filter program to process and filter an allelic data-set.
	mode = 1;
	fill_pformat_data2_command(command, NUMPROCS, bmap, bped, phe, phenoName, 
							   NUMPROCS, run_id, LDIR, RDIR, rowMajor, nperm, mode, mind, geno, maf,
							   person_missingness, founders, hwe);

	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// 2nd call to join the fragments and create the filtered data-set.
	// This call is always a single core operation.
	mode = 2;
	fill_pformat_data2_command(command, 1, bmap, bped, phe, phenoName, 
							   NUMPROCS, run_id, LDIR, RDIR, rowMajor, nperm, mode, mind, geno, maf,
							   person_missingness, founders, hwe);

	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// First call to 'pepi_qt_v2.exe', no permuted calculation spread across N-cores.
	mode = MODE_NON_PERMUTED;
	fill_pepi_qt_v2_command(command, NUMPROCS, run_id, samples, filtered_bmap, filtered_sped, LDIR, RDIR, nperm, mode);
	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// 2nd call on a single core is to combine the results from all ranks and apply a p-value threshold to find
	// putative interesting interactions in the data-set.
	mode = MODE_NON_PERMUTED_MERGE;
	fill_pepi_qt_v2_command(command, 1, run_id, samples, filtered_bmap, filtered_sped, LDIR, RDIR, nperm, mode);
	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// 3rd call on N-cores to run the permutation.
	// each core runs a randomised data-set so running epistasis on large data-set is
	// not advised as it will take ages.
	mode = MODE_PERMUTED;
	fill_pepi_qt_v2_command(command, NUMPROCS, run_id, samples, filtered_bmap, filtered_sped, LDIR, RDIR, nperm, mode);
	printf("\n");
	result = system(command);
	printf("\n");
	CU_ASSERT(result == 0);
	if (result != 0) return;

	// 4th call on single core to merge the permutation results into a single file.
	// Single core operation
	mode = MODE_PERMUTED_MERGE;
	fill_pepi_qt_v2_command(command, 1, run_id, samples, filtered_bmap, filtered_sped, LDIR, RDIR, nperm, mode);
	result = system(command);
	CU_ASSERT(result == 0);
	if (result != 0) return;

	printf("\nFile Outputs:-\n");
	printf("\nrun_ex_0016.fepi.csv = Interesting Interactions (QT) with permutation statistics\n");
	printf("run_ex_0016.fslog.txt = log file containing filter settings\n");
	printf("run_ex_0016.rseed.txt = Initialisation value used by PRNG\n");
	echo_spacer(); 
}
